1 and 2. Look at the chapter on interactive graphics and, specifically, the code to display a subject's MRICloud data as a sunburst plot. Do the following. Display this subject's data as a Sankey diagram. Display as many levels as you can for type = 1, starting from the intracranial volume. Put this in a file called hw4.ipynb.

import pandas as pd
import plotly.express as px
import numpy as np
import os
import plotly.graph_objects as go

## load in the hierarchy information
url = "https://raw.githubusercontent.com/bcaffo/MRIcloudT1volumetrics/master/inst/extdata/multilevel_lookup_table.txt"
multilevel_lookup = pd.read_csv(url, sep = "\t").drop(['Level5'], axis = 1)
multilevel_lookup = multilevel_lookup.rename(columns = {
    "modify"   : "roi", 
    "modify.1" : "level4",
    "modify.2" : "level3", 
    "modify.3" : "level2",
    "modify.4" : "level1"})
multilevel_lookup = multilevel_lookup[['roi', 'level4', 'level3', 'level2', 'level1']]
multilevel_lookup.head()

	roi	level4	level3	level2	level1
0	SFG_L	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L
1	SFG_R	SFG_R	Frontal_R	CerebralCortex_R	Telencephalon_R
2	SFG_PFC_L	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L
3	SFG_PFC_R	SFG_R	Frontal_R	CerebralCortex_R	Telencephalon_R
4	SFG_pole_L	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L

## Now load in the subject data
id = 127
subjectData = pd.read_csv("kirby21AllLevels.csv")
subjectData = subjectData.loc[(subjectData.type == 1) & (subjectData.level == 5) & (subjectData.id == id)]
# Show only this subject's region (ROI) and volume 
subjectData = subjectData[['roi', 'volume']]
## Merge the subject data with the multilevel data
subjectData = pd.merge(subjectData, multilevel_lookup, on = "roi")
subjectData = subjectData.assign(icv = "ICV")
subjectData = subjectData.assign(comp = subjectData.volume / np.sum(subjectData.volume))
subjectData.head()

	roi	volume	level4	level3	level2	level1	icv	comp
0	SFG_L	12926	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L	ICV	0.009350
1	SFG_R	10050	SFG_R	Frontal_R	CerebralCortex_R	Telencephalon_R	ICV	0.007270
2	SFG_PFC_L	12783	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L	ICV	0.009247
3	SFG_PFC_R	11507	SFG_R	Frontal_R	CerebralCortex_R	Telencephalon_R	ICV	0.008324
4	SFG_pole_L	3078	SFG_L	Frontal_L	CerebralCortex_L	Telencephalon_L	ICV	0.002227

subjectData = subjectData.drop(['volume'], axis=1)

df2=subjectData.groupby(['icv', 'level1'], as_index=False)['comp'].count()
df2.head()
df3=subjectData.groupby(['level1', 'level2'], as_index=False)['comp'].count()
df3.head()
df4=subjectData.groupby(['level2','level3'], as_index=False)['comp'].count()
df4.head()
df5=subjectData.groupby(['level3','level4'], as_index=False)['comp'].count()
df5.head()
df6=subjectData.groupby(['level4','roi'], as_index=False)['comp'].count()
df6.head()

	level4	roi	comp
0	AG_L	AG_L	1
1	AG_R	AG_R	1
2	ALIC_L	ALIC_L	1
3	ALIC_R	ALIC_R	1
4	Amyg_L	Amyg_L	1

df2.columns = ['a', 'b', 'comp']
df3.columns = ['a', 'b', 'comp']
df4.columns = ['a', 'b', 'comp']
df5.columns = ['a', 'b', 'comp']
df6.columns = ['a', 'b', 'comp']
df7 = df2.append(df3)
df8 = df7.append(df4)
df9 = df8.append(df5)
df10 = df9.append(df6)
df10

	a	b	comp
0	ICV	CSF	32
1	ICV	Diencephalon_L	8
2	ICV	Diencephalon_R	8
3	ICV	Mesencephalon	8
4	ICV	Metencephalon	20
...	...	...	...
269	midbrain_R	Midbrain_R	1
270	midbrain_R	RedNc_R	1
271	midbrain_R	Snigra_R	1
272	post_DPWM_L	PCR_L	1
273	post_DPWM_R	PCR_R	1

500 rows × 3 columns

import chart_studio.plotly as py

def genSankey(df,cat_cols=[],value_cols='',title='Sankey Diagram'):
    # maximum of 6 value cols -> 6 colors
    colorPalette = ['#FFD43B','#646464','#4B8BBE','#306998']
    labelList = []
    colorNumList = []
    for catCol in cat_cols:
        labelListTemp =  list(set(df[catCol].values))
        colorNumList.append(len(labelListTemp))
        labelList = labelList + labelListTemp
        
    # remove duplicates from labelList
    labelList = list(dict.fromkeys(labelList))
    
    # define colors based on number of levels
    colorList = []
    for idx, colorNum in enumerate(colorNumList):
        colorList = colorList + [colorPalette[idx]]*colorNum
        
    # transform df into a source-target pair
    for i in range(len(cat_cols)-1):
        if i==0:
            sourceTargetDf = df[[cat_cols[i],cat_cols[i+1],value_cols]]
            sourceTargetDf.columns = ['source','target','count']
        else:
            tempDf = df[[cat_cols[i],cat_cols[i+1],value_cols]]
            tempDf.columns = ['source','target','count']
            sourceTargetDf = pd.concat([sourceTargetDf,tempDf])
        sourceTargetDf = sourceTargetDf.groupby(['source','target']).agg({'count':'sum'}).reset_index()
        
    # add index for source-target pair
    sourceTargetDf['sourceID'] = sourceTargetDf['source'].apply(lambda x: labelList.index(x))
    sourceTargetDf['targetID'] = sourceTargetDf['target'].apply(lambda x: labelList.index(x))
    
    # creating the sankey diagram
    data = dict(
        type='sankey',
        node = dict(
          pad = 15,
          thickness = 20,
          line = dict(
            color = "black",
            width = 0.5
          ),
          label = labelList,
          color = colorList
        ),
        link = dict(
          source = sourceTargetDf['sourceID'],
          target = sourceTargetDf['targetID'],
          value = sourceTargetDf['count']
        )
      )
    
    layout =  dict(
        title = title,
        font = dict(
          size = 10
        )
    )
       
    fig = dict(data=[data], layout=layout)
    return fig

import plotly

fig = go.Figure(genSankey(df10,
                cat_cols=['a','b'],
                value_cols='comp',
                title='Sankey Diagram of Brain Regions'),
               )
fig.update_layout(
    autosize=False,
    width=1000,
    height=2000,
    paper_bgcolor="LightSteelBlue",
)
#plotly.offline.plot(fig, validate=False)
fig.show()

#fig.to_html()

1. Link to website: https://eyuasa.github.io/

1. Create the opioid sqlite database from https://smart-stats.github.io/ds4bio_book/book/_build/html/sqlite.html. However, only go to the step where the csv files are read into the database. Then exit sqlite and you should have a file opioid.db that has the data. Next, read the three tables into pandas dataframes and do the data wrangling from the sqlite chapter directly in pandas. Add the python code to your hw4.ipynb file.

import sqlite3 as sq3
import pandas as pd

con = sq3.connect("opioid.db")
# cursor() creates an object that can execute functions in the sqlite cursor

sql = con.cursor()

population = pd.read_sql_query("SELECT * from population", con)
land = pd.read_sql_query("SELECT * from land", con)
annual = pd.read_sql_query("SELECT * from annual", con)

# you have to close the connection
con.close

<bound method NDFrame.head of                BUYER_COUNTY BUYER_STATE countyfips STATE COUNTY  \
0          1        AUTAUGA          AL      01001     1      1   
1          2        BALDWIN          AL      01003     1      3   
2          3        BARBOUR          AL      01005     1      5   
3          4           BIBB          AL      01007     1      7   
4          5         BLOUNT          AL      01009     1      9   
...      ...            ...         ...        ...   ...    ...   
28260  28261       WASHAKIE          WY      56043    56     43   
28261  28262         WESTON          WY      56045    56     45   
28262  28263        SKAGWAY          AK      02230     2    230   
28263  28264  HOONAH ANGOON          AK      02105     2    105   
28264  28265     PETERSBURG          AK      02195     2    195   

         county_name                               NAME    variable  year  \
0            Autauga            Autauga County, Alabama  B01003_001  2006   
1            Baldwin            Baldwin County, Alabama  B01003_001  2006   
2            Barbour            Barbour County, Alabama  B01003_001  2006   
3               Bibb               Bibb County, Alabama  B01003_001  2006   
4             Blount             Blount County, Alabama  B01003_001  2006   
...              ...                                ...         ...   ...   
28260       Washakie           Washakie County, Wyoming  B01003_001  2014   
28261         Weston             Weston County, Wyoming  B01003_001  2014   
28262        Skagway       Skagway Municipality, Alaska  B01003_001  2014   
28263  Hoonah Angoon  Hoonah-Angoon Census Area, Alaska  B01003_001  2014   
28264     Petersburg         Petersburg Borough, Alaska  B01003_001  2014   

      population  
0          51328  
1         168121  
2          27861  
3          22099  
4          55485  
...          ...  
28260       8444  
28261       7135  
28262        996  
28263       2126  
28264       3212  

[28265 rows x 11 columns]>

land.head()

	Areaname	STCOU	LND010190F	LND010190D	LND010190N1	LND010190N2	LND010200F	LND010200D	LND010200N1	LND010200N2	...	LND110210N1	LND110210N2	LND210190F	LND210190D	LND210190N1	LND210190N2	LND210200F	LND210200D	LND210200N1	LND210200N2
0	UNITED STATES	00000	0	3787425.08	0000	0000	0	3794083.06	0000	0000	...	0000	0000	0	251083.35	0000	0000	0	256644.62	0000	0000
1	ALABAMA	01000	0	52422.94	0000	0000	0	52419.02	0000	0000	...	0000	0000	0	1672.71	0000	0000	0	1675.01	0000	0000
2	Autauga, AL	01001	0	604.49	0000	0000	0	604.45	0000	0000	...	0000	0000	0	8.48	0000	0000	0	8.48	0000	0000
3	Baldwin, AL	01003	0	2027.08	0000	0000	0	2026.93	0000	0000	...	0000	0000	0	430.55	0000	0000	0	430.58	0000	0000
4	Barbour, AL	01005	0	904.59	0000	0000	0	904.52	0000	0000	...	0000	0000	0	19.59	0000	0000	0	19.61	0000	0000

5 rows × 34 columns

annual.head()

	BUYER_COUNTY	BUYER_STATE	year	count	DOSAGE_UNIT	countyfips
0	ABBEVILLE	SC	2006	877	363620	45001
1	ABBEVILLE	SC	2007	908	402940	45001
2	ABBEVILLE	SC	2008	871	424590	45001
3	ABBEVILLE	SC	2009	930	467230	45001
4	ABBEVILLE	SC	2010	1197	539280	45001

population.head()

		BUYER_COUNTY	BUYER_STATE	countyfips	STATE	COUNTY	county_name	NAME	variable	year	population
0	1	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2006	51328
1	2	BALDWIN	AL	01003	1	3	Baldwin	Baldwin County, Alabama	B01003_001	2006	168121
2	3	BARBOUR	AL	01005	1	5	Barbour	Barbour County, Alabama	B01003_001	2006	27861
3	4	BIBB	AL	01007	1	7	Bibb	Bibb County, Alabama	B01003_001	2006	22099
4	5	BLOUNT	AL	01009	1	9	Blount	Blount County, Alabama	B01003_001	2006	55485

annual = pd.read_csv("county_annual.csv") population = pd.read_csv("county_pop_arcos.csv") land = pd.read_csv("land_area.csv")

annual.describe()

		BUYER_COUNTY	BUYER_STATE	year	count	DOSAGE_UNIT	countyfips
count	27758	27758	27758	27758	27758	27758	27758
unique	27758	1864	57	9	12212	24890	3040
top	14165	WASHINGTON	TX	2006	1	100	NA
freq	1	273	2032	3100	50	32	760

population.describe()

		BUYER_COUNTY	BUYER_STATE	countyfips	STATE	COUNTY	county_name	NAME	variable	year	population
count	28265	28265	28265	28265	28265	28265	28265	28265	28265	28265	28265
unique	28265	1849	51	3142	51	325	1850	3140	2	9	24017
top	14165	WASHINGTON	TX	47137	48	3	Washington	NA	B01003_001	2006	7588
freq	1	279	2286	9	2286	441	279	15	28250	3142	6

land.describe()

		Areaname	STCOU	LND010190F	LND010190D	LND010190N1	LND010190N2	LND010200F	LND010200D	LND010200N1	...	LND110210N1	LND110210N2	LND210190F	LND210190D	LND210190N1	LND210190N2	LND210200F	LND210200D	LND210200N1	LND210200N2
count	3198	3198	3198	3198	3198	3198	3198	3198	3198	3198	...	3198	3198	3198	3198	3198	3198	3198	3198	3198	3198
unique	3198	3193	3198	4	3152	4	1	4	3147	3	...	1	1	4	1943	4	1	4	1949	3	1
top	2721	Franklin, VA	18035	0	0	0000	0000	0	0	0000	...	0000	0000	0	0	0000	0000	0	0	0000	0000
freq	1	2	1	3190	5	3191	3198	3191	5	3194	...	3198	3198	3190	26	3191	3198	3191	26	3194	3198

4 rows × 35 columns

population = population.drop(columns="", axis=1)

annual = annual.iloc[: , 1:]
#annual = population.drop(columns="", axis=1)
#land = population.drop(columns="", axis=1)

land = land.iloc[: , 1:]

population.head()

		BUYER_COUNTY	BUYER_STATE	countyfips	STATE	COUNTY	county_name	NAME	variable	year	population
0	1	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2006	51328
1	2	BALDWIN	AL	01003	1	3	Baldwin	Baldwin County, Alabama	B01003_001	2006	168121
2	3	BARBOUR	AL	01005	1	5	Barbour	Barbour County, Alabama	B01003_001	2006	27861
3	4	BIBB	AL	01007	1	7	Bibb	Bibb County, Alabama	B01003_001	2006	22099
4	5	BLOUNT	AL	01009	1	9	Blount	Blount County, Alabama	B01003_001	2006	55485

print(annual[annual['countyfips'] == "NA"])

      BUYER_COUNTY BUYER_STATE  year count DOSAGE_UNIT countyfips
187       ADJUNTAS          PR  2006   147      102800         NA
188       ADJUNTAS          PR  2007   153      104800         NA
189       ADJUNTAS          PR  2008   153       45400         NA
190       ADJUNTAS          PR  2009   184       54200         NA
191       ADJUNTAS          PR  2010   190       56200         NA
...            ...         ...   ...   ...         ...        ...
27753           NA          NV  2007   447      200600         NA
27754           NA          NV  2008     5        2200         NA
27755           NA          OH  2006    23        5100         NA
27756           NA          PR  2006    10       17800         NA
27757           NA          PR  2007     2        1300         NA

[760 rows x 6 columns]

print(annual[(annual['countyfips']== "NA") & (annual["BUYER_STATE"] != "PR")])

      BUYER_COUNTY BUYER_STATE  year count DOSAGE_UNIT countyfips
10071         GUAM          GU  2006   319      265348         NA
10072         GUAM          GU  2007   330      275600         NA
10073         GUAM          GU  2008   313      286900         NA
10074         GUAM          GU  2009   390      355300         NA
10075         GUAM          GU  2010   510      413800         NA
...            ...         ...   ...   ...         ...        ...
27751           NA          MA  2006   247      114900         NA
27752           NA          NV  2006   380      173600         NA
27753           NA          NV  2007   447      200600         NA
27754           NA          NV  2008     5        2200         NA
27755           NA          OH  2006    23        5100         NA

[74 rows x 6 columns]

annual["countyfips"] = annual["countyfips"].astype(str)

annual.loc[(annual["BUYER_STATE"] == "AR") & (annual["BUYER_COUNTY"] == "MONTGOMERY"), "countyfips"] = "05097"

print(annual.loc[(annual["BUYER_STATE"] == "AR") & (annual["BUYER_COUNTY"] == "MONTGOMERY")])

      BUYER_COUNTY BUYER_STATE  year count DOSAGE_UNIT countyfips
17429   MONTGOMERY          AR  2006   469      175390      05097
17430   MONTGOMERY          AR  2007   597      241270      05097
17431   MONTGOMERY          AR  2008   561      251760      05097
17432   MONTGOMERY          AR  2009   554      244160      05097
17433   MONTGOMERY          AR  2010   449      247990      05097
17434   MONTGOMERY          AR  2011   560      313800      05097
17435   MONTGOMERY          AR  2012   696      339520      05097
17436   MONTGOMERY          AR  2013   703      382300      05097
17437   MONTGOMERY          AR  2014   491      396900      05097

#sqlite> delete from annual where BUYER_COUNTY = "NA"
print(annual[annual['BUYER_COUNTY'] == "NA"])

      BUYER_COUNTY BUYER_STATE  year count DOSAGE_UNIT countyfips
27741           NA          AE  2006     2         330         NA
27742           NA          CA  2006    47       12600         NA
27743           NA          CT  2006   305       78700         NA
27744           NA          CT  2007   112       30900         NA
27745           NA          CT  2008    48       15000         NA
27746           NA          FL  2006     9         900         NA
27747           NA          FL  2007     7         700         NA
27748           NA          GA  2006   114       51700         NA
27749           NA          IA  2006     7        2300         NA
27750           NA          IN  2006   292       39300         NA
27751           NA          MA  2006   247      114900         NA
27752           NA          NV  2006   380      173600         NA
27753           NA          NV  2007   447      200600         NA
27754           NA          NV  2008     5        2200         NA
27755           NA          OH  2006    23        5100         NA
27756           NA          PR  2006    10       17800         NA
27757           NA          PR  2007     2        1300         NA

annual = annual[annual['BUYER_COUNTY'] != "NA"]
print(annual[annual['BUYER_COUNTY'] == "NA"])

Empty DataFrame
Columns: [BUYER_COUNTY, BUYER_STATE, year, count, DOSAGE_UNIT, countyfips]
Index: []

# sqlite> create table land_area as select Areaname, STCOU, LND110210D from land;
land_area = land[["Areaname", "STCOU", "LND110210D"]]

# sqlite> alter table land_area rename column STCOU to countyfips;
land_area = land_area.rename(columns={"STCOU": "countyfips"})
land_area

	Areaname	countyfips	LND110210D
0	UNITED STATES	00000	3531905.43
1	ALABAMA	01000	50645.33
2	Autauga, AL	01001	594.44
3	Baldwin, AL	01003	1589.78
4	Barbour, AL	01005	884.88
...	...	...	...
3193	Sweetwater, WY	56037	10426.65
3194	Teton, WY	56039	3995.38
3195	Uinta, WY	56041	2081.26
3196	Washakie, WY	56043	2238.55
3197	Weston, WY	56045	2398.09

3198 rows × 3 columns

#sqlite> create table county_info as select * from population left join land_area using(countyfips);
county_info = pd.merge(population, land_area, left_on = "countyfips", right_on = "countyfips")

county_info.head()

	BUYER_COUNTY	BUYER_STATE	countyfips	STATE	COUNTY	county_name	NAME	variable	year	population	Areaname	LND110210D
0	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2006	51328	Autauga, AL	594.44
1	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2007	52405	Autauga, AL	594.44
2	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2008	53277	Autauga, AL	594.44
3	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2009	49584	Autauga, AL	594.44
4	AUTAUGA	AL	01001	1	1	Autauga	Autauga County, Alabama	B01003_001	2010	53155	Autauga, AL	594.44

len(land)

3198

len(land_area)

3198

len(county_info)

28256

len(population)

28265

1. Create an interactive scatter plot of average number of opiod pills by year plot using plotly. See the example here. Don't do the intervals (little vertical lines), only the points. Add your plot to an html file with your repo for your Sanky diagram and host it publicly. Put a link to your hosted file in a markdown cell of your hw4.ipynb file. Note, an easy way to create a webpage with this graphic is to export an ipynb as an html file.

annual.head()

	BUYER_COUNTY	BUYER_STATE	year	count	DOSAGE_UNIT	countyfips
0	ABBEVILLE	SC	2006	877	363620	45001
1	ABBEVILLE	SC	2007	908	402940	45001
2	ABBEVILLE	SC	2008	871	424590	45001
3	ABBEVILLE	SC	2009	930	467230	45001
4	ABBEVILLE	SC	2010	1197	539280	45001

annual['year'] = pd.to_numeric(annual['year'])

<ipython-input-65-b23c1c0c5aea>:1: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
  annual['year'] = pd.to_numeric(annual['year'])

annual['DOSAGE_UNIT'] = pd.to_numeric(annual['DOSAGE_UNIT'])

<ipython-input-66-2b4db4e52f15>:1: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
  annual['DOSAGE_UNIT'] = pd.to_numeric(annual['DOSAGE_UNIT'])

annual_group = annual.groupby(by=['year'])['DOSAGE_UNIT'].mean().reset_index()

annual_group

	year	DOSAGE_UNIT
0	2006	2.654669e+06
1	2007	2.995906e+06
2	2008	3.254271e+06
3	2009	3.526038e+06
4	2010	3.783656e+06
5	2011	4.035583e+06
6	2012	3.993064e+06
7	2013	3.861752e+06
8	2014	3.768738e+06

fig_opioids = px.scatter(annual_group, x="year", y="DOSAGE_UNIT", title="Average Opioid Pill Shipments by Year")
fig_opioids.show()

fig_opioids_2 = px.line(annual_group, x="year", y="DOSAGE_UNIT", title="Average Opioid Pill Shipments by Year", markers=True)
fig_opioids_2.show()

Link to website: https://eyuasa.github.io/

import sqlite3 as sq3 import pandas as pd

conn = sq3.connect("opioid.db")

sql = con.cursor()

annual.to_sql(str(annual), conn, if_exists='append', index=False)

population.to_sql(str(population), conn, if_exists='append', index=False) land.to_sql(str(land), conn, if_exists='append', index=False)